A Call Option equals:
\begin{align*} C(K) = e^{-r(T-t)}\mathbb{E}(S_T - K|S_t)^+ = S_t \Phi(d_1) - Ke^{-r(T - t)}\Phi(d_2) \end{align*}A Put Option equals:
\begin{align*} P(K) = e^{-r(T-t)}\mathbb{E}(K-S_T|S_t)^+ = Ke^{-r(T - t)}\Phi(-d_2) - S_t \Phi(-d_1) \end{align*}$\Phi(x)$ is the standard normal cumulative density function. The probability that $X \le .3$ is $\mathrm{P}({X \le .3}) = \Phi(.3) = 0.62 $
import pandas as pd
import numpy as np
import scipy.stats as stats
import plotly.express as px
import plotly.graph_objects as go
import plotly.io as pio
import sympy
pio.renderers.default = 'notebook'
class Contract:
def __init__(self,S0,K,R,Sigma,T):
"""
Class defines GBM contract
S0: Spot price at time 0
K: Strike price
R: risk-free interest rate
Sigma: annualized volatility of the option
T: Time to maturity
"""
self.S0 = S0
self.K = K
self.R = R
self.Sigma = Sigma
self.T = T
class CALL(Contract):
def __init__(self,S0,K,R,Sigma,T):
Contract.__init__(self,S0,K,R,Sigma,T)
class PUT(Contract):
def __init__(self,S0,K,R,Sigma,T):
Contract.__init__(self,S0,K,R,Sigma,T)
class Option_Pricer:
def __init__(self,contract):
"""
contract: use the CALL or PUT derived classes to defined the contract
"""
self.contract = contract
self.S0, self.K, self.R, self.Sigma, self.T = contract.S0,contract.K,contract.R,contract.Sigma, contract.T
self.Type = self.contract.__class__.__name__
def d1(self):
D1 = (np.log(self.S0/self.K)+(self.R+ .5*self.Sigma**2)*self.T)/(self.Sigma*np.sqrt(self.T))
return D1
def d2(self):
d_1 = self.d1()
D2 = d_1 - self.Sigma*np.sqrt(self.T)
return D2
def pv(self):
return self.K*np.exp(-self.R*self.T)
def __str__(self):
return (f"{self.Type} option contract with strike {self.K} and volatility {self.Sigma}.")
def price_derivative(self):
D1 = self.d1()
D2 = self.d2()
PV = self.pv()
if self.Type == 'CALL':
return self.S0*stats.norm.cdf(D1)- PV*stats.norm.cdf(D2)
else:
return PV*stats.norm.cdf(-D2)-self.S0*stats.norm.cdf(-D1)
def plot_contract(self):
prices = self.price_derivative()
if len(prices) <= 1:
raise Exception('Cannot Graph 1 Price.')
DF = pd.DataFrame(prices, index = self.S0, columns = [f"{self.Type}(K = {self.K})"])
DF.index.name = 'Spot Price'
fig = px.line(DF).update_layout(yaxis_title = 'Derivative Price')
return fig
def plot_delta(self):
title = f"{self.Type}(K = {self.K}, T = {self.T},vol = {self.Sigma},r = {self.R})"
Deltas = self.delta()
DF = pd.DataFrame(Deltas,index = self.S0, columns = [f"{self.Type}(K = {self.K})"])
DF.index.name = 'Spot Price'
fig = px.line(DF,title = title).update_layout(yaxis_title = "Delta")
fig.add_shape(type ='line',x0 = self.K,x1 = self.K,y0= np.where(self.Type == 'CALL',0,-100),y1=np.where(self.Type == "CALL",100,0), line = dict(color = "red"),name = 'Strike')
return fig
def plot_gamma(self):
title = f"{self.Type}(K = {self.K}, T = {self.T},vol = {self.Sigma},r = {self.R})"
gammas = self.gamma()
DF = pd.DataFrame(gammas,index = self.S0, columns = [f"{self.Type}(K = {self.K})"])
DF.index.name = 'Spot Price'
gamma_max = max(gammas)
fig = px.line(DF,title = title).update_layout(yaxis_title = "Gamma")
fig.add_shape(type ='line',x0 = self.K,x1 = self.K,y0=0,y1=gamma_max+1, line = dict(color = "red"),name = 'Strike')
return fig
def plot_vega(self):
title = f"{self.Type}(K = {self.K}, T = {self.T},vol = {self.Sigma},r = {self.R})"
vegas = self.vega()
DF = pd.DataFrame(vegas,index = self.S0, columns = [f"{self.Type}(K = {self.K})"])
DF.index.name = 'Spot Price'
vega_max = max(vegas)
fig = px.line(DF,title = title).update_layout(yaxis_title = "Vega")
fig = fig.add_shape(type ='line',x0 = self.K,x1 = self.K,y0=0,y1=vega_max+1, line = dict(color = "red"),name = 'Strike')
return fig
def delta(self):
if self.Type == 'CALL':
return stats.norm.cdf(self.d1())*100
else:
return -stats.norm.cdf(-self.d1())*100
def gamma(self):
return stats.norm.pdf(self.d1())/(self.S0*self.Sigma*np.sqrt(self.T))*100
def vega(self):
return self.S0*np.sqrt(self.T)*stats.norm.pdf(self.d1())
def return_values(self):
price = self.price_derivative()
delta =self.delta()
gamma = self.gamma()
vega = self.vega()
vals = np.array([price,delta,gamma,vega])
df = pd.DataFrame([vals],columns = ['Price','Delta','Gamma','Vega'],index = [f'{self.Type}(K = {self.K})'])
return df
S0 = np.linspace(80,120)
K = 100
T = .5
r = .02
sigma = .3
call1 = Option_Pricer(CALL(S0,K,r,sigma,T))
put1 = Option_Pricer(PUT(S0,K,r,sigma,T))
$\Delta = \frac{\partial V}{\partial S_{t}}$
If a put contract has a $\Delta = -.5$ this means for every unit increase in the price of the underlying the option price decreases by $-.5 \cdot (\Delta_{s}$). If $S$ increases by $\$5$, then the price of the put option theoretically declines by $\$2.5$ or $\$250$ notional.
Below are two plots that show how the delta of a call and put option changes with respect to their spot prices.
call1.plot_delta()
put1.plot_delta()
Gamma is the second partial derivative of any derivative contract with respect to its price. It essentially measures the change in the option's delta.
Delta hedging works well for small price movements and for options with greater time to maturity, however, as the markets become more dislocated and maturity decreases, gamma can drastically impact profit or loss.
The gamma of a long call or put is always positive, where as the gamma of a short call or short put is always negative.
If a call option has a $\Gamma = .03$ and $\Delta = .55$ and the underlying asset increases by $\$3$, then the theoretical price increase in the call option equals $0.55 \cdot 3 + \frac{0.03 \cdot 3^{2}}{2} = 1.65 + 0.13 = 1.785$.
call1.plot_gamma()
If $\upsilon = 12.3$, then a $1\%$ increase in the implied volatility increases the option price by $0.01 \cdot 12.3 = .123 \ \text{or by} \ \$12.3 \ \text{Notional}$.
\begin{align*} \boxed{\text{Put Vega} = \text{Call Vega} = \upsilon = \frac{\partial C}{\partial \sigma} = S \sqrt T N'(d_1)} \end{align*}The vega for a long call or put option is always positive.
put1.plot_vega()
call1.plot_vega()
Let's say we believe implied volatility is trading lower so we buy a straddle on the SPY ETF Index. Let's say we want to buy a 0-Delta, $\$1000$ vega notional option strategy. This is what a straddle accomplishes because we are hedged against price movements but long implied volatility.
There are two-options being sold right now:
Call(S0 = 430,K =426,R = .03,$\sigma = .13$,$\tau = .25$) = $\$15.0159 \cdot 100$.
Put(S0 = 430,K =433,R = .03,$\sigma = .13$,$\tau = .25$) = $\$11.028 \cdot 100$
| Price | $\Delta$ | $\Gamma$ | $\upsilon$ | |
|---|---|---|---|---|
| CALL(K = 426) | 15.016 | 0.614729 | 0.013679 | 82.2008 |
| PUT(K = 433) | 11.0281 | -0.483679 | 0.0142615 | 85.7008 |
What combination of the Call and the Put option contracts will delta hedge my position and profit one-thousand dollars when implied volatility increases by one percent? Find the number of puts and calls to buy such that our portfolio will have the these dynamics: \begin{align*} \Delta_{p} &= 0 \\ \upsilon_{p} &= 1000 \end{align*}
If the implied volatility increase by $1\%$, then my option portfolio vega increases by $0.01 \cdot 1000 = \$10$ which is $\$1000$ notional value.
This system of equations is equal to the matrix equation:
S1 = 430
K_c = 426
r = .03
vol = .13
Tau = .25
K_p = 433
call_2 = Option_Pricer(CALL(S1,K_c,r,vol,Tau))
put_2 = Option_Pricer(PUT(S1,K_p,r,vol,Tau))
DF = pd.concat([call_2.return_values(),put_2.return_values()],axis =0)
DF
| Price | Delta | Gamma | Vega | |
|---|---|---|---|---|
| CALL(K = 426) | 15.015958 | 61.472942 | 1.367904 | 82.200768 |
| PUT(K = 433) | 11.028142 | -48.367873 | 1.426148 | 85.700800 |
Solving the problem we get:
| CALL(K = 426) | PUT(K = 433) | |
|---|---|---|
| Delta | 61.4729 | -48.3679 |
| Vega | 82.2008 | 85.7008 |
| Position | 5.23197 | 6.65032 |
We need to buy about 5.23 calls and buy 6.65 puts. Seldomly, you'll ever get a perfect whole number for a solution, so just round to the nearest one!
NEW = (DF.T).loc[['Delta','Vega'],:]
A = sympy.Matrix([[61.47,-48.36],[82.2,85.7]])
P_g = sympy.Matrix([0,1000])
wts = A.solve(P_g)
NEW.loc['Position'] = [float(wts[0]),float(wts[1])]
A = sympy.Matrix([[61.47,-48.36],[82.2,85.7]])
P_g = sympy.Matrix([0,1000])
wts = A.solve(P_g)